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I. INTRODUCTION 



The basic reasons for using an autopilot on a ship are 
to reduce the number of operator controls , minimize rudder 
orders for Course-keeping , and reduce the traveling time 
between destinations. 

In the past autopilots were used extensively to maintain 
a given course. They work well but most autopilots are for 
course control only, the effect of an ocean current and/or 
wind during travel are not considered. When a disturbance is 
applied, the ship follows the desired course but does not 
follow the given track without correction. This problem has 
been studied [Ref. 1]. The correction needed to compensate 
for errors due to the effect of a current and/or wind were 
found . 

The purpose of this thesis is to develop a procedure to 
compensate for the error in position and reduce it as much 
as possible, then the autopilot can be used for both 
Course-keeping and Track following. 

Today, many ships have computers on board. The computer 
can be used to solve the problem of Track following. 
First, the desired trajectories are stored in the computer. 
When the computer receives measured X and Y positions the 
computer will provide the desired position with respect to X 
and Y and compare with the actual position. An error signal 
is obtained from trajectory inf ormation , and can be used to 
drive the ship close to the desired trajectory. 

Actual position must be provided by the navigation 
system. The navigation system used in the Track following 
autopilot must be accurate , must give continuous information 
about the position of the ship , and the system must be 
useable everywhere in the world. It is felt that the 
NAV STAR/ GLOBAL POSITIONING SYSTEM can be used. 
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The NAVSTAR/ GLOBAL POSITIONING SYSTEM ( GPS ), currently 
being developed by the Department of Defense, is a satellite 
based navigation system that will provide the user with 
extremely accurate three-dimensional position, velocity and 
time information on a 24-hour basis and in all weather 
conditions at any point on the earth. The user position is 
determined by measuring its range to four satellites . For 
more details refer to [Refs. 2,3], 



10 



II. COURSE-KEEPING 



A. THE PRIMARY TRANSFER FUNCTION OF THE SYSTEM 

The equations of motion by Davidson [Ref. 4 ] are: 

+ c^ip - mS 2 = c^<S 

nft + c Q - c ib =* c <$ 

K m y 

where Q(s)= (Z/V)d 

0 = turning angular velocity ; 

V = ship speed; 
ip = drift angle ; 

6 = helm angle ; 

Z = ship length ; 

denotes d/ds = U/V)d/dt , s = (V/Ot 
m = (mi-c £ ) 

mi,m2,n = coefficient of inertia ; 

c c c c = coefficient of resistance ; 

Z’. k* m’ f 

c c = coefficient of rudder force. 

P* X 

The former equation is the equation of lateral transla- 
tion and the latter is the equation of turning angular 
motion. 

Then rewriting the equations of motion in terms of time 
and taking the Laplace transform, we obtain the turning rate 
of the ship in steering [Ref. 5 ] is: 

0 ( s ) _ K( I+T3 s ) 6 ( s ) + fTiT?s+(Ti+T?)T 9(0 -)+TiT26(0-) 

"( 1 +Tis) ( 1 +T 2 s) ( 1 +Tis) ( 1 +T 2 s) 

The first term corresponds to the ship motion excited by 
the steering, and the second corresponds to the memory of her 
motion at the beginning of the steering. Therefore, a rela- 
tional function is: 
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( 1 ) 



0*( s ) = K(l+T,s) 
S(s) (1+Tis) (1+T 2 s) 



Equation(l) describes the response character of the ship 
to steering, which may be called the transfer function of 
the ship in steering. 

In this study of automatic Course-keeping we will be 
working with a 200,000 DWT super- tanker of the following 



characteristics [Ref. 6], 




length = 310.00 


meters 


breadth = 47.16 


meters 


Draft = 18.90 


meters 


Steering Quality indices 


: . 


T1 = -269.3 


seconds 


T2 = 9.3 


seconds 


T3 = 20.0 


seconds 



K = -0.0434 red/sec 

Maximum Rudder Deflection = 30 degree 

Maximum Rudder rate = 2.32 degrees / second 

Substitute Tl,T2,T3 and K in equation (1), we get: 

0(s) _ -0. 0434(l+20s) 

6 ( s ) (1-269. 3s) (1+9. 3si 



( 2 ) 



B. STABILIZING THE SYSTEM 

From equation(2) ,we see that one pole of the steering 
transfer function is in the right half plane, so the system 
is unstable. We have to stabilize the system by using an 
autopilot . 
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Assuming the autopilot has a gain G, Figure 2.1 is the 
block diagram of the system with an autopilot. 




Figure 2.1 Block Diagram Of The System With An Autopilot. 



The open loop transfer function for the system is 

H(s) GK( 1 + T a s ) 

" s(1+Tis)(1+T 2 s) (1+T e s) 

Where — = the rudder servo transfer function. 

<1+T e s) 

From experience the value of T E is 1 to 2 seconds and a 
good choice is 1.7 seconds, so the only parameter available 
to make the system stable is G. 

To get the values of G for stability, we refer to the 
Routh criterion: 



GK( 1+T3s ) . -1 

s(l+T e s)(l+Tls)(l+T2s) 

s(l + Te s) (1 + Tls) (1 + T2s) +GK( 1 + T3s ) = 0 

s(l+1.7s)(l-269.3s)(l+9.3s)-0.0434G(l+20s) = 0 

s‘‘ + 0.692s 3 +0.06s 2 -(2.35Xl0* 4 -2.04xl0~‘‘G)s+1.02X10‘ 5 G = 0 

s 4 1 0.06 1.02X10' 5 G 
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s 3 0.692 - (2. 35X1 0* 4 - 2.04* 10^* G) 0 

s 2 A 1 . 02'xlO* 5 G 

s 1 B 0 

s° 1 . 02X10* 5 G 



Where 

A = 0.692 0 . 06 + (2 . 35x10*** - 2 . 04x1 0‘ 5 G ) 

0.692 

B = -A(2.35Xl0~ 4 -2.04xl0~ 4 G)-0. 692 1.02xl0~ 6 G 

A 



To find the limits of G for stability , the values of A,B 
and 1.02xl0' 5 G must be greater than zero. 

In the s 2 row 
A > 0 

0.692 0.06+(2.35Xl0* 4 -2.04X10* 4 G) > 0 

4 . 1755X10' 2 -2 . 04X10" 4 G > 0 
G < 204.68 



In the s 1 row 
B > 0 

-A(2.35X10' 4 -2.04xl0' 4 G)-0. 692X1. 02X10' 5 G > 0 
After manipulation. 

G 2 - 85G + 240 < 0 
3 < G < 82 
In the s 1 row 

1 . 02X10* 5 G > 0 

G > 0 

So the condition of G for stability is 
3 < G < 82 

By using DSL/360, the system was simulated. The system 
was represented by the block diagram of Figure 2.1. The 
computer program for this system is contained in Appendix A. 
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Values of G were selected between 12 and 36, which are the 
ones -with better time constants. 

Figure 2.2 to Figure 2.9 are the computer outputs for 
different values of G. By comparing, we found that when 
G=24.2 (Figure 2. 4, 2. 5) and G=30 (Figure 2. 6, 2. 7) the 
settling time is almost the same and shorter, with less 
oscillation than for other values of G. When we compare the 
rudder angle for both values, we can see that when G=30 the 
rudder angle is larger so in this case G = 24.. 2 is the best 
value . 

Next we simulate a disturbance that produces the rate of 
turn caused by waves. For a large super tanker the rate of 
turn can easily reach as much as 0.2 or 0.3 
deg ./ sec . (0 . 00349 or 0.00524 rad. /sec.). 

In this thesis the rate of turn 0.2 deg. /sec. was used 
in the program. Figure 2.10 to Figure 2.17 are the computer 
outputs of this program. Figure 2.12,2.13 show that when 
G=24.2 the settling time and rudder angle are better than 
other values. Figure 2.18 is the Bode plot of this system 
which phase margin = 11.17 degrees and gain margin =-18.46 
DB that are lower than normally used in practise. 
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Figure 2 . 2 



Ship Heading(l) and Heading Command (2) 
with G=20. 
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Figure 2.3 



Rudder An^le(l^ and Ship Heading(2) 
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Figure 2.4 Ship Heading(l) and Heading Command(2) 

with G=24.2. 
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Figure 2.5 



Rudder Anglefl) and Ship Heading(2) 
with G=24 . 2 . 
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Figure 2.6 Ship Heading(l) and Heading Command(2) 

with G=30. 
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Figure 2.7 



Rudder Angle(l) and Ship Heading(2) 
with G=30. 
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Figure 2.8 Ship Heading(l) and Heading Command(2) 

with G=36 . 
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Figure 2.9 Rudder Angle(l) and Ship Heading(2) 

with G= 36 . 



23 



cP.OO 4.00 8.00 12.00 16.00 20.00 24.00 23.00 32.00 36.00 



HERD AND HERDRF VS TIME 
HE ADRF= 1 0 IN2=0. 00349 G=20.0 




X SCALE = 100.00 UNITS/INCH RUN NO. I 

Y SC RL E = 4.00 UNITS/INCH PLOT NO. I 



Figure 2 . 10 



Ship Heading(l) and Heading Command(2) 
with Disturbance and G=20. 
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Figure 2.11 Rudder Angle(l) and Ship Heading(2) 
with Disturbance and G=20. 
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Figure 2.12 _ Ship Heading ( 1 ) and Heading Command(2) 

with Disturbance and G=24.2. 
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Figure 2.13 Rudder Angle(l) and Ship Heading(2) 
with Disturbance and G=24.2. 
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Figure 2.15 Rudder Angle(l) and Ship Heading(2) 
with Disturbance and G=30. 
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Figure 2.16 



Ship Heading(l) and Heading Command(2) 
with Disturbance and G=36. 
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Figure 2.17 Rudder Angle(l) and Ship Heading(2) 
with Disturbance and G=36. 
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Figure 2.18 Bode Plot of the System with Steering Servo. 
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C. IMPROVEMENT OF THE SYSTEM 



Now the system is stable but the rudder angle is too 
large which is not practical at all. It is necessary to 

improve the system performance by using a compensator to 
obtain acceptable transient and rudder angle performances. 

1 . Cascade Lead Compensation 

Lead filter compensation has a transfer function as 

follows : 

1 (s + Z) _ P . 2l • _ 1 + ^ s 

a (s+P) Z P 

Where oc=Z/P and P>Z. The system is shown in Figure 

2.19. 



(l+-p-s) l+-ps 



0I1TWO.QAWC t 
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e 



i 
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Figure 2.19 The System with the Compensator. 

The characteristic equation is: 

(G/ot)K(l + T3s) (s + z) 

+ s ( s+p ) (1 + Tls) (1 + T2s) (1 + Tg s) 

-0 . 0434(G/<x ) (s+ p) ( l + 20s ) 

' + s(s+p)(l-269.3s)(l + 9.3s)(l + 1.7s) = ° 
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0 



s (s + p )( 1- 269. 3s ) (1+9. 3s) (1+1. 7s) -0.0434 (G/a) (s+a P ) (l+20s) 

4257. 633s 5 +(2946. 49+4257. 63 3p)s 4 +(258. 3+2946. 49p)s 3 + 

{- 1+258. 3p+.868(G/(X)}s 2 + {(-l+.868G)p+ .0434 (G/a)} s+.0434Gp 

= 0 

From this equation, the NPS PAROLE program was used 
to find the best values of P and Z. The family of root loci 
are given on Figure 2.20. 

From experience a good value of (X was 0.1 and we 
select the damping ratio about 0.42. We obtain: 

P = 0.4 
ot = 0.1 
Z = CXP = 0.04 

Hence the transfer function of the lead filter 

A 

compensator is: 



1 ' f (T7^4 s 

1+ TiV 



1 + 2 5s 
1+2. bs 



Figure 2.21-2.24 are the computer outputs of Figure 
2.19 with different values of G. It can be observed that to 
have a quick response a high gain is necessary but this 

could demand an excessive operation of the rudder which can 

be bad for the following reasons: 

(1) To get a quick response , faster rudder operation and 
greater rudder angle are needed. 

(2) This faster response may cause accelerations that 
are too sharp for personnel aboard ship. 

(3) Too much rudder action would cause the rudder 

machine with all the hardware to have its life 

reduced . 

In Figure 2.22 - 2.24 the system reaches steady 

state in very short time but the rudder angle is too high. 
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REAL AXIS IUN I TS PER INCH) = 0.2500 
IMAG AXIS (UNITS PER INCH) = 0.2500 

ALPHA=. 001 , . 01 , . 1 
VIMUKTANANDA G=24.2 




Figure 2.20 Family of the Root Loci. 
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In Figure 2.21 the rudder angle is at the limit but the 
system takes a long time to reach steady state. 

So, it was necessary to reduce the rudder angle and 
simultaneously reduce the settling time of the system. In 
order to introduce a limiting value of rudder rate, the 
steering gear transfer function must be replaced in the 
system. The block diagram for an equivalent circuit was 
given on Figure 2.25. 

Where K3 =1/T e = 0.588 

With this transformation , we have the system on 
Figure 2.26 . 

Figure 2.27 - 2.30 are the computer outputs of this 
block diagram with the gain selected before. We can see that 
when the value of the gain(Kl) is greater than 10 , the system 
reaches steady state in a short time (about 3 minutes for 
Kl=10) but the rudder angle reaches the limit on one side. 
In Figure 2 . 30 ,K1=3 . 5 , the rudder angle is within the limits 
but the settling time is still too long. Then we must try to 
find a value of gain(Kl) between 3.5 and 10 to get the best 
results for rudder angle and settling time. 

Figure 2.31 - 2.33 are the outputs of this trial 
with gain: 4. 6, 5. 3 and 6. We can see that: 

For gain = 6 the settling time is about 450 
sec . ( 7 . 5min . ) the maximum rudder rudder angle is about 28 
deg. and the heading error about 0.77 deg.. 

For gain = 5.3 the settling time is about 470 
sec . ( 7 . 8min . ) the maximum rudder angle is about 26 deg. and 
the heading error about 0.9 deg.. 

For gain = 4.6 the settling time is about 
510sec . (8 . 5min) the maximum rudder angle is about 22 deg. 
and the heading error about 1.1 deg.. 
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Figure 2.21 Rudder Angle(l) and Ship Heading(2) 
without Limiter and G=3.5. 
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Figure 2.22 Rudder Angle(l) and Ship Heading(2) 
without Limiter and G=10. 
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Figure 2.23 Rudder Angle (1) and Ship Heading (2) 
without Limiter and G = 15. 
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Figure 2.24 Rudder Anele(l) and Ship Heading(2) 
without Limiter and G=24.2T & 1 
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Figure 2.25 



Block Diagram of the Steering Gear. 
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Figure 2.26 The Complete Block Diagram for Course-Keeping. 

So we can conclude that the best value for 
Course-keeping with the limiter is 4.6. Figure 2.34 is the 
Bode plot of this system with phase margin = 42.2 degrees 
and gain margin = -15.97 DB that are better than the Bode 
plot of Figure 2.18. Then the system of Course-keeping is 
satisfactory . 
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Figure 2.27 Rudder Angle(l) and Ship Heading(2) 
with Limiter and Kl=24.2. 
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Figure 2.28 Rudder Angle(l) and Ship Heading(2) 
with Limiter and Kl=15. 
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Figure 2.29 Rudder Angle(l) and Ship Heading(2) 
with Limiter and Kl=10. 
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Figure 2.30 Rudder Angle(l) and Ship Heading(2) 
with Limiter and Kl=3.5. 



45 



,- 25.00 - 20.00 - 15.00 - 10.00 - 5.00 0.00 5.00 10.00 15.00 20.00 



RUDDER RNGLE 4 HERD VS TIME 
I N2 = 0 . 00349 K 1 =4 . 6 HERDRF = 3 




U.00 



80.00 



160.00 



240.00 



320.00 



400.00 



480.00 



560.00 



640.00 



XSCRLE= 80.00 UNITS/INCH 

YSCflLE= 5.00 UNITS/INCH 



RUN NO. 1 
PLOT NO. 1 



Figure 2.31 Rudder Angle(l) and Ship Heading(2) 
with Limiter ana Kl=4.6. 
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Figure 2.32 Rudder Angle(l) and Ship Heading(2) 
with Limiter and Kl=5.3. 
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Figure 2.33 Rudder Angle(l) and Ship Heading(2) 
with Limiter and Kl=6. 
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Figure 2.34 Bode Plot of the System 
with Steering Servo and Filter. 
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III. TRACK FOLLOWING 



Using the ship heading controller obtained previ- 
ously, and assuming that the velocity of the ship is 
constant (at 14 knots) during underway on the open sea, the 
position of the ship can be found at any time. 

In order to find the position of the ship, a right hand 
rectangular coordinate system is established, the origin of 
which is chosen to be in the body itself, as shown in Figure 
3.1. 




Figure 3.1 Orientation of the Space Axis(X0,Y0) 
and the Moving Axis(X,Y). 

The origin and the axes are fixed with respect to the 
body but movable with respect to other systems of coordinate 



50 



axes fixed in space.lt is assumed that the two systems coin- 
cide at t=0. 

The transformation from the ship to space coordinate 
system is defined by the following relations , obtained from 
Figure 3 . 1 







X = VcosG + UsinQ 






Y* = Vs in© - UcosQ 


and 




X = X0+ /xdt 






Y = Y0+/Ydt 


where 


X 


= Velocity in X-direction 




Y 


= Velocity in Y-direction 




XO 


= Initial position of X 




YO 


= Initial position of Y 




V 


= Ship velocity 




U 


= Lateral velocity 


Assuming 


constant velocity and no 


equations 


become : 



lateral force, the 



X = VcosQ 
Y = Vsin9 



and 

X = dx/dt 
Y = dy/dt 



then 

X = XO + J dx 
Y = YO + Jdy 

Knowing initial values of X and Y,the coordinates of the 
ship are calculated. Figure 3.2 is the block diagram of this 
procedure . 
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Figure 3.2 Course-Keeping and Coordinate Calculation. 

These equations were included in the program of the 
Course-keeping autopilot in order to find the position of 
the ship. 

Figure 3. 3-3. 4 are the output of this program. 

The purpose of the Track following is to keep the ship 
following the desired trajectory from the beginning to the 
destination. When the ship is not on the track a course 
correction must be calculated to return the ship to the 
track. The corrected course is a function of the distance to 
the track. 

A. FINDING THE COURSE CORRECTION 

Initially the navigator must design the trajectory from 
the beginning point to the destination on the map. He has to 
know the course and speed. When the ship is underway, it is 
necessary to know(measure) the position of the ship. 
Although the course of the ship may be the same as the 
desired course, it may not be on the desired trajectory due 
to the effect of sea current and wind. Therefore the navi- 
gator has to calculate a new course to keep the ship 
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Ship Heading(l) and Desired Heading(2). 
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Figure 3.4 Distance in X(l) and Y(2). 
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following the desired trajectory until the ship reaches the 
destination. This situation is depicted in Figure 3.5. 




Figure 3.5 The Trajectory Followed by the Ship. 

The procedure of an autopilot for track-following is the 
same as that used by a navigator. First the coordinates of 
the desired trajectory are stored in the computer as shown 
in Figure 3.6 The position of the ship is measured by 
NAVSTAR/GPS . The computer calculates the errors in X and Y 
positions and calculates the course correction. These errors 
should be zero in order to keep the ship always on the 
desired trajectory. 

The procedure to calculate the course correction(9 c ) is 
shown in Figure 3.7 
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Figure 3.6 The Coordinate of The Desired Trajectory. 




Figure 3.7 The Algorithm to Compute The Course Correction. 

DS/Xe = SIN 9 

DS = Xe SIN 9 

m = Vds 1 + £ 

SIN 6c = DS/M 

0c = SIN 1 (DS/M) 
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B. FINDING THE DESIRED COURSE 



The desired trajectory is stored in the computer, the 
desired course can be found in a simple way, Figure 3.8 shows 
the method. The X axis is divided into constant small inter- 
vals. For each point of X,the corresponding Y coordinate is 
found. The desired course is obtained by taking the arctan- 
gent of the result of the division between increment in Y 
and increment in X. 




Figure 3.8 Algorithm to Compute The Desired Course. 

0 , - 

©2 = TAN/ Y 2 -Y| \ 

vx^-x,; 

03 = TAnV y 3 -Y 2 \ 

VX3-X2J 

04 ~ TAN*/ Y 4 -Y 3 \ 

VX4-X3J 

05 = TAnVy 5 -Y4\ 

XX5-X4/ 
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For selecting the value of the constant small interval 
of X,we consider Figure 3.9 




Figure 3.9 The Condition such that the Algorithm Fails. 

If the value of X is big, the straight line AB is not on 
the curve, so the desired course is not true. To get the true 
desired heading the straight line AB has to be on the curve 
all the time. To accomplish this we have to select the value 
of X small enough, so the straight line AB is part of the 
curve. 

C. FIRST ATTEMPT TO DO THE TRACK FOLLOWING 

From the Figure 3.2, X and Y positions are obtained. 
With these values, the desired heading and trajectory are 
obtained. By using an integrator with the gain K4(as shown 
in the block diagram of figure 3.2 ) we will obtain a new 
block diagram which is shown in Figure 3.10. The purpose of 
the integrator is to eliminate steady state position errors. 

A suitable value of integrator gain is K4=.01.The value 
was found by repeated simulation runs. Figure 3.11-3.14 are 
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Figure 3.10 Block Diagram of An Autopilot 
With The Integrator. 

the results obtained without disturbances. From Figure 
3.11-3.12 the actual trajectory and the desired trajectory 
are coincident. In Figure 3.14 the time for the output to 
reach steady state is about 800 seconds(13.3 minutes). 

The system was also simulated with disturbances. Figure 
3.15-3.18 are the results with K4=.01. In Figure 3.15-3.16 
the actual trajectory and the desired trajectory are not 
coincident but the distance between them gradually increases 
with time. In Figure 3.18 the time for the output to reach 
steady state is about 1600 seconds (26 . 67 minutes). 
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Figure 3.11 Desired Trajectory and X when K4=.01. 



60 



isimi 




DESIRED TRAJECTORY AND 
ACTUAL TRAJECTORY ( Y I VS TIME 



54 

46 



4? 



56 



31 



24 



18 



12 * 



8 



I 



Figure 3.12 Desired Trajectory and Y when K4=.01. 
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Figure 3.13 Rudder Angle when K4=.01. 
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SHIP HERDING AND DESIREO HEflOING VS TlflE 
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Figure 3.14 Ship Heading and Heading Command when K4=.01. 
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Figure 3.15 Desired Trajectory and X with Disturbance 

when K4 = . 01 . 
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Figure 3.16 Desired Trajectory and Y with Disturbance 

when K4= . 01 . 



65 



isanw A 



RUOHNG 1 DEGREES 1 




RUOOER ANGLE VS TIME 



Figure 3.17 Rudder Angle with Disturbance 

when K4= . 01 . 
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Figure 3.18 Ship Heading and Heading Command 
with Disturbance when K4 =.01. 
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Additional simulation runs included disturbances , but it 
was also assumed that the ship was not on the desired 
trajectory at the beginning. Figure 3.19-3.22 are the 
results. In Figure 3. 19-3. 20, we can notice that the actual 
trajectory stays away from the desired trajectory all the 
time, there are no control actions to put the ship back to 
desired trajectory. 

The next step was to include the course correction(0 ) 
in the control loop, as shown in Figure 3.23. 

Simulations were run with course correction and gain 
K5 , a best value of K5=.75 was obtained. Figure 3.24-3.27 
are the results using this value. In Figure 3.24-3.25 the 
ship comes back to the desired trajectory and follows it all 
the time. The time from the initial position to the desired 
trjectory is about 1400 seconds(23.3 minutes). 
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Figure 3.19 Desired Trajectory and X with 
Initial Value of X=l. 
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Figure 3.20 Desired Trajectory and Y with 
Initial Value of x=l. 
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Figure 3.21 Rudder Angle with Initial Value of x=l. 
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SHIP HERDING flNO DESIRED HERO I MG VS TIME 
IN2-0. 00349 X-1 K4-.01 



Figure 3.22 Ship Heading and Heading Command with 

Initial Value of x= I. 
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Figure 3.23 The Complete Block Diagram for Course-Keeping 

and Track-Following. 
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DESIRED TRAJECTORY AND 
ACTUAL TRAJECTORY (Y) VS TIME 



Figure 3.24 Desired Trajectory and 
Actual Traj ectory (y ) with K5=.75. 
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Figure 3.25 Desired Trajectory and 
Actual Trajectory (x) with K5=.75. 
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Figure 3.26 Rudder Angle with K5=.75. 
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Figure 3.27 Ship Heading and Desired Heading 

with K= . 7 5 . 
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IV. CONCLUSION 



Use of the computer to control the autopilot for both 
Course-keeping and Track following achieves an accuracy for 
steady state conditions as high as may be desired. An 
optimal value of gain for the autopilot was obtained by 
using the PAROLE program to study families of Root Loci. 
Both negative and positive feedback can be used, but only 
negative feedback is recommended. 

In Chapter 2 for Course-keeping , the Bode plot is used to 
show that the filter POLE ZERO and gain were the best for 
the system. The actual heading in steady state was exactly 
the desired heading with zero error. 

In Chapter 3 the problem of acquiring the track was 
studied. If the Track following mode is turned on when the 
ship is many ship lengths off track, the amount of rudder 
activity is determined by the distance of f- track, the nature 
of the Tracking algorithm , and the shape of track itself, as 
well as the initial value of the ship heading. It was found 
that the initial ship heading was important , and rudder 
activity could be minimized by proper choice of the initial 
heading angle. This suggests that further study of the 
effect of initial heading angle is desirable , since manual 
adjustment of ship heading prior to activating the autopilot 
could greatly reduce rudder motions. Once the ship reached 
the track, the system was able to follow the desired trajec- 
tory exactly in calm and rough sea. 

Further studies may be conducted to include the cost 
function into these programs to minimize the fuel consump- 
tion and to apply optimal control theory to find the optimal 
track for ship maneuvering at sea. 
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APPENDIX A 



COMPUTER PROGRAM 



* THIS PROGRAM REFERS TO FIGURES 2.2-2.17 
//SOMM07 JOB ( 2259 , 1435 ) , * PROJECT 1 , CLASS=B 
//*MAIN ORG=NPGVMl . 2259P 

//* FORMAT PR, DDNAME=PLOTX. SYSVECTR, DEST=LOCAL 
// EXEC DSL 
//DSL. INPUT DD * 

TITLE AUTOPILOT DESIGN WITHOUT COMPENSATION 
INTGER NPLOT 

CONST NPLOT=l , IC1=0 . , IC2=0 . , IC3=0 . , IC4=0 . , Cl=3 . 1415927 



PARAM HEADRF 
PARAM IN2 
PARAM G 
PARAM K 
PARAM TE 
PARAM T1 
PARAM T2 
PARAM T3 
DERIVATIVE 
INI 
ERROR 
DELTAR 
DELTA 
AMPLI 
COMP 
CORREC 
THETA 1 
THETA 
HEAD 
RUDANG 
ERRORD 



10 . 

0.0 

20 

-.0434 

1.7 

-269.3 

9.3 

20 . 

HEADRF*C1/180. 

INI-THETA 

G*ERROR 

REALPL( IC1,TE, DELTAR) 
K*DELTA 

LED LAG ( IC2,T3,T1, AMPLI) 
REALPL ( IC3 , T2 , COMP ) 
IN2+CORREC 
INTGRL (14, THETA1 ) 

THETA* 180. /Cl 
DELTA* 180. /Cl 
HEADRF- HEAD 



SAMPLE 

CALL DRWG(1, 1, TIME, HEAD) 

CALL DRWG( 1,2, TIME, HEADRF) 
CALL DRWG(2, 1, TIME, RUDANG) 
CALL DRWG( 2 , 2 , TIME , HEAD) 
TERMINAL 

CALL ENDRW( NPLOT) 

PRINT 8. , HEADRF, HEAD, ERRORD, RUDANG 
CONTRL FINTIM=300. , DELT= . 8 , DELS=8 . 
END 
STOP 

//PLOT. PLOTP ARM DD * 

&PLOT SCALE= . 65 SEND 
//PLOT. SYS IN DD * 

HEAD AND HEADRF VS TIME 
HEADRF=10 IN2=0 . 0 G=20 
RUDDER ANGLE AND HEAD VS TIME - 
HEADRF=10 IN2=0 . 0 G=20 
/* 

/* 
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FILE: EE 44181 DSL 



A1 



* THIS PROGRAM REFERS TO FIGURES 2.21-2.24 
//SOMMAR3 JOB ( 2259 , 1435 ), 1 EE44181 CLASS=B 
//*MAIN ORG=NPGVMl . 2259P , LINES=( 10 ) 

//* FORMAT PR / DDNAME=PLOTX. SYSVECTR , DEST= LOCAL 
// EXEC DSL 
//DSL. INPUT DD * 

TITLE AUTOPILOT DESIGN WITH COMPENSATION 
INTGER NPLOT 

CONST NPL0T=4 , IC1=0. , 102=0., IC3=0. , IC4=0 . ,01=3.1415927 
PARAM HEADRF=0 . 

PARAM IN2=0. 00349 
PARAM G=3 . 5 
PARAM K=- . 0434 
PARAM TE=1 . 7 



PARAM Tl=-269 . 3 
PARAM T2=9 „ 3 
PARAM T3=20 . 



DERIVATIVE 

INI 

ERROR 

FILTER 

DELTAR 

DELTA 

AMPLI 

COMP 

CORREC 

THETA1 

THETA 

HEAD 

RUDANG 

ERRORD 



HEADRF*C1/18Q . 

INI-THETA 

LED LAG ( 0 . ,25. ,2. 5, ERROR) 
G* FILTER 

REALPL( IC1,TE, DELTAR) 
K*DELTA 

LED LAG ( IC2 , T3 , T1 , AMPLI ) 
REALPL( IC3 , T2 , COMP ) 
IN2+C0RREC 
INTGRL (14, THETA1 ) 

THETA* 180. /Cl 
DELTA* 180. /Cl 
HEADRF-HEAD 



SAMPLE 



CALL DRWG( 1,1, TIME, RUDANG) 
CALL DRWG (1,2, TIME , HEAD ) 
TERMINAL 

CALL ENDRW( NPLOT) 

PRINT 3 ., HEAD , RUDANG 

CONTRL FINTIM=800 . , DELT= . 8 , DELS=8 . 

END 

PARAM G=10 
END 

PARAM G=15 
END 

PARAM G=24 . 2 

END 

STOP 

//PLOT. PLOTP ARM DD * 

&PLOT SCALE= . 65 SEND 
//PLOT. SYS IN DD * 

RUDDER ANGLE AND HEAD VS TIME 
SOMMART G=3 . 5 IN2=. 00349 
RUDDER ANGLE AND HEAD VS TIME 
SOMMART G=10 IN2=. 00349 
RUDDER ANGLE AND HEAD VS TIME 
SOMMART G=15 IN2=. 00349 
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FILE: EE44185 DSL 



A1 



* THIS PROGRAM REFERS TO FIGURES 2.27-2.33 
//SOMMAR7 JOB (2259 , 1435 ),' EE44184 ' , CLASS=B 
//♦MAIN ORG=NPGVMl . 2259P , LINES=( 10 ) 

//♦FORMAT PR, DDNAME=PLOTX. SYSVECTR, DEST= LOCAL 
// EXEC DSL 
//DSL. INPUT DD * 

TITLE AUTOPILOT DESIGN WITH COMPENSATION AND LIMITER 
PARAM IN2=0 . 00349 
PARAM Kl=24. 2 

PARAM K2=-. 0434, K3=. 588, .. . 

Pl=- . 1222 , P2= . 1222 , P3=- . 52 , P4= . 52 , K4=0 . 0 
PARAM HEADRF=3 
INTGER NPLOT 

CONST NPL0T=1,C1=3. 1415927 
DERIVATIVE 

INI = HEADRF*C1/180 . 

ERROR1 = INI -THETA 

FILTER = LEDLAG(0 . , 25 . ,2.5, ERROR1) 

DELREF = FILTER*K1 

DELDES = LIMIT(P3,P4, DELREF) 

DELTA1 = DELDES -DELTAC 
DELDOT = LIMIT (PI, P2,DELTA1) 

DELDOC = K3*DELDOT 
DELTAC = INTGRL ( 0 . , DELDOC ) 

AMPLI = K2*DELTAC 

COMP = LEDLAG(0. ,20. , -269 . 3 , AMPLI ) 

CORREC = REALPL ( 0 . ,9.3, COMP ) 

THETA 1 = IN2+CORREC 
THETA = INTGRL ( 0 . , THETA 1 ) 

RUDANG = DELTAC* 180. /Cl 
HEAD = THETA* 180. /Cl 
ERRORD = HEADRF-HEAD 

SAMPLE 

CALL DRWG(1, 1, TIME, RUDANG) 

CALL DRWG( 1,2, TIME, HEAD) 

TERMINAL 

CALL ENDRW( NPLOT) 

PRINT 8 . , RUDANG, HEAD , ERRORD 
CONTRL FINTIM=700 . , DELT= . 7 , DELS=7 . 

INTEG RKSFX 

END 

STOP 

//PLOT . PLOTPARM DD * 

&PLOT SCALE= . 65 &END 
//PLOT . SYSIN DD * 

RUDDER ANGLE & HEAD VS TIME 
IN2=0. 00349 Kl=24. 2 HEADRF=3 
/* 

/* 
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FILE: T2 



DSL 



A1 



* THIS PROGRAM REFERS TO FIGURES 3.3-3. 4 
//SOMMA21 JOB ( 2259 , 1435 ) , 1 PROJECT ' , CLASS=B 
//*MAIN ORG=NPGVM1.2259P,LINES=(10) 

//* FORMAT PR , DDNAME=PLOTX . SYSVECTR , DEST=LOCAL 
// EXEC DSL 
//DSL. INPUT DD * 

TITLE SIMULATION OF SHIP DYNAMICS (TEST OF FOLLOWING) 
PARAM K1=4.6,K2=-.0434,K3=.588, . . . 

Pl=- . 1222 , P2= . 1222 , P3=- . 52 , P4= . 52 
PARAM V=0. 0038889 

* SPEED IN MILES/SECOND WHICH IS EQUIVALENT TO 14 KNOTS 
PARAM IN2=0 . 0 , HEADRF=5 

CONST NPLOT=l , Cl=3 . 1415927 

INTGER NPLOT 

DERIVATIVE 

INI = HEADRF*C1/180 . 

El = INI -THETA 

FILTER= LED LAG ( 0 . , 25 . , 2 . 5 , El ) 

E2 = FILTER*K1 

E3 = E1*K4 

E4 = INTGRL ( 0 . , E3 ) 

DELREF= E2+E4 

DELDES= LIMIT(P3 / P4 / DELREF) 

DELTA1= DELDES-DELTAC 
DELDOT= LIMIT ( PI , P2 r DELTA1 ) 

DELDOC= K3*DELDOT 
DELTAC= INTGRL ( 0 . # DELDOC ) 

AMPLI = K2*DELTAC 

COMP = LEDLAG(0. ,20. / -269.3 / AMPLI) 

CORREC= REALPL ( 0 . , 9 . 3 , COMP ) 

THETA 1= IN2+CORREC 
THETA = INTGRL ( 0 . , THETA 1 ) 

HEAD = THETA *180. /Cl 
RUDANG= DELTAC* 180 . /Cl 
DX = V*COS( THETA) 

DY = V* SIN (THETA) 

X = INTGRL ( 0. ,DX) 

Y = INTGRL ( 0 . , DY) 

SAMPLE 

CALL DRWG( 1,1, TIME, HEAD) 

CALL DRWG( 1,2, TIME, HEAD RF) 

CALL DRWG ( 2 , 1 , T I ME , X ) 

CALL DRWG (2, 2, TIME, Y) 

TERMINAL 

CALL ENDRW( NPLOT) 

PRINT 8 . , HEAD , HE ADRF , X , Y 

CONTRL FINTIM=650. , DELT= . 4, DELS= . 65 

INTEG RKSFX 

END 

STOP 

//PLOT . PLOTP ARM DD * 

&PLOT SCALE= . 65 SEND 
//PLOT. SYS IN DD * 

HEAD AND HE ADRF VS TIME 
HEADRF=5 IN2=0 
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FILE: S0M2 



DSL 



A1 



* THIS PROGRAM REFERS TO FIGURES 3.11-3.22 

TITLE SIMULATION OF SHIP DYNAMICS (TEST OF FOLLOWING) 
PARAM K1=4.6,K2=-. 0434, K3=. 588 , K4=0 . 01, . . . 

Pl=- . 1222 , P2=. 1222 , P3=- . 52 , P4= . 52 
PARAM V=0. 0038889 

* SPEED IN MILES/SECOND WHICH IS EQUIVALENT TO 14 KNOTS 
PARAM EPSILX=0. 125,0=-. 074 

PARAM IN2=0 . 00349 
CONST Cl=3. 1415927 

NLFGEN YPATH=-20. / -30.,-10.,-12.5 / 0.,0. , 10 . , 10 . , 20 . , . . . 

18. . 30. .25. .40. .31. .50. .36. . 60. .40. .70. .43 . .80. . . . . 

45.. 90.. 46. 5. 100.. 48. 

NLFGEN XPATH=-30 . ,-20. ,-12.5,-10. ,0. ,0. ,10. ,10. ,18. , . . . 

20 . . 25 . . 30 . . 31 . . 40 . . 36 . . 50 . . 40 . . 60 . . 43 . . 70 . . 45 . . . . . 

80.. 46. 5. 90.. 48.. 100. 

DYNAMIC 

TRAJE1 = NLFGEN ( YP ATH , X ) 

TRAJE2 = NLFGEN (XP ATH, Y) 

DELTAX = (X+EPSILX) -X . 

DELTAY = NLFGEN (YP ATH, X+EP SI LX) -NLFGEN ( YP ATH, X) 

Z = DELTAY/DELTAX . . _ 

ZR = ATAN(Z) 

THECOM = ZR* 180 . /Cl ' 

DERIVATIVE 

INI = THECOM*C1/180. 

E = INI-THETA 

El = 10. *E 

FILTER = ZEROPL(0. , .04, 0.4, El) 

* E2 = FILTER*K1 

E3 = E*K4 

E4 = INTGRL ( 0 . , E3 ) 

DELREF = E2+E4 

DELDES = LIMIT(P3,P4, DELREF) 

DELTA1 = DELDES -DELTAC 
DELDOT = LIMIT ( PI ,P2,DELTA1) 

DELDOC = K3*DELDOT 
DELTAC = I NTGRL ( 0. , DELDOC) 

AMPLI1 = K2*DELTAC 
AMPLI = G*AMPLI 1 

COMP = ZEROPL(0. , .05, -.0037, AMPLI) 

CORREC = REALPL ( 0 . ,9.3, COMP ) 

THETA 1 = IN2+CORREC 

THETA2 = INTGRL(0. 7853982, THETA1) 

PROCED THETA = COR ( THETA2 , Cl ) 

THE T A3 = THETA2 

3 IF(ABS(THETA3) .LT.C1) GO TO 5 
I F ( THETAS . LT . 0 ) GO TO 4 
THETA3 = THETA2-C1 

GO TO 3 

4 THETA3 = THETA2+C1 
GO TO 3 

5 THETA = THETA2 
ENDPRO 

HEAD = THETA *180. /Cl % 

RUDANG = DELTAC *180. /Cl 
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FILE : S0M2 



DSL 



A1 



DX 

DY 

X 

Y 



= V*COS( THETA) 

= V* SIN (THETA) 

= INTGRL(0 . ,DX) 
= INTGRL(0 . , DY) 



PRINT 125. ,RUDANG,HEAD,THEC0M,X,TRAJE1,Y,TRAJE2 
CONTRL FINTIM=24000. , DELT=6 . , DELS=20 . 

SAVE ( G1 ) 12 . , HEAD , THECOM 
SAVE ( G2 ) 12 . , TRAJE1 , Y 
SAVE (G3 ) 12 . , TRAJE2 , X 
SAVE ( G4 ) 12 . , RUDANG 

GRAPH ( Gl/Gl , DE=TEK618 , P0=0 , . 5 ) TIME ( LE=8 . 0 , SC=500 . , . . . 

NI=10,UN=’SEC ' ), . . . 

HEAD ( LE=9 ,NI=9 , LO=37, SC=2 . , UN= ’ DEGREES ’ ) , . . . 

THECOM (LE=9,NI =9, PO=8. , LO=37,SC=2. , UN= ’ DEGREES ' ) 
GRAPH (G2/G2,DE=TEK618,PO=0, .5) TIME ( LE=8 . , SC=2000 . , . . . 

NI = 12 , UN= 1 SEC' ) , . . . 

TRAJE1 ( LE=9 , NI=9 , L0=0 , SC=6 . , UN=' MILES ' ) , . . . 

Y ( LE=9 , NI=9 , PO=8 . / LO=0,SC=6. ,UN=' MILES' ) 

GRAPH (G3/G3 , DE=TEK618 , P0=0 , . 5 ) TIME ( LE=8 . , SC=2000 . ... . 

NI=12 ,UN=' SEC' ) , . . . 

TRAJE2 (LE=9 y NI=9 , L0=0 , SC=9 . , UN=' MILES ' ) ... . 

X(LE=9,NI=9,PO=8. ,LO=0,SC=9. ,UN=' MILES' ) 

GRAPH ( G4/G4, DE=TEK618 ) TIME ( LE=8 . , SC=300 . , . . . 

NI=12,UN='SEC' ), . . . 

RUDANG( LE=9 ,NI=9 , L0=-30 / SC=7 . y UN= 'DEGREES ' ) 

LABEL (Gl) SHIP HEADING AND DESIRED HEADING VS TIME 

LABEL (Gl) IN2=0. 00349 X=0 K4=.01 

LABEL (G2 ) DESIRED TRAJECTORY AND 

LABEL (G2 ) ACTUAL TRAJECTORY (Y) VS TIME 

LABEL (G3 ) DESIRED TRAJECTORY AND 

LABEL ( G3 ) ACTUAL TRAJECTORY (X) VS TIME 

LABEL (G4) RUDDER ANGLE VS TIME 

END 

STOP 
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FILE: S0M1 



DSL 



A1 



* THIS PROGRAM REFERS TO FIGURE 3.24 TO 3.27 

TITLE SIMULATION OF SHIP DYNAMICS (TEST OF FOLLOWING) 
PARAM Kl=4 . 6 , K2=- . 0434, K3=. 588 , K4=0 . 01 , K5=0 . 75 , . . . 

Pl=- . 1222 , P2= . 1222 , P3=- . 52 , P4=. 52 
PARAM V=0. 003 8889 

* SPEED IN MILES/SECOND WHICH IS EQUIVALENT TO 14 KNOTS 
PARAM EPSILX=0 . 125 , G=- . 074 

PARAM L=2 . 5 , IN2=0 . 00349 
CONST Cl=3. 1415927 

NLFGEN YPATH=-20 . , -30 . , -10 . , -12 . 5 , 0 . , Q. , 10 . , 10 . , 20 . , . . . 

18. . 30. .25. .40. .31. .50. .36. .60. .40. .70. .43. ... . 

80 . . 45 . . 90 . .46.5. 100 . . 48 . 

NLFGEN XPATH=-30. , -20. , -12 . 5 , -10 . , 0. , 0 . , 10 . , 10 . , 18. , . . . 

20. . 25. .30. .31. .40. .3. 6. ,50. ,40. ,60. ,43. ,70. ,.. . 

45 . , 80 . , 46 . 5 , 90 . , 48 . , 100 . 

DYNAMIC 

TRAJE1 = NLFGEN ( YPATH , X ) 

TRAJE2 = NLFGEN ( XP ATH , Y ) 

DELTAX = (X+EPSILX) -X 

DELTAY = NLFGEN ( YPATH ,X+EP SI LX) -NLFGEN ( YPATH ,X) 

2 = DELTAY/DELTAX 

2R = AT AN ( Z ) 

THECOM = ZR*180./C1 
DERIVATIVE 

INI = THECOM*C1/180 . 

E = IN1-THETA+THETAC 

El = 10. *E 

FILTER = ZER0PL(0. , ..04,0. 4, El) 

E2 = FILTER*K1 “ 

E3 = E*K4 

E4 = INTGRL(0. , E3 ) 

DELREF = E2+E4 

DELDES = LIMIT (P3 , P4, DELREF) 

DELTA1 = DELDES-DELTAC 
DELDOT = LIMIT ( PI , P2 , DELTA1 ) 

DELDOC = K3*DELDOT 
DELTAC = INTGRL(0. , DELDOC) 

AMPLI1 = K2*DELTAC 
AMPLI = G*AMPLI1 

COMP = ZEROPL(0. , .05, -.0037, AMPLI) 

CORREC = REALPL ( 0 . ,9.3, COMP ) 

THETA 1 = IN2+CORREC 

THETA2 = INTGRL( 0 . 9599311 , THETA1 ) 

PROCED THETA = COR ( THETA2 , Cl ) 

THETA3 = THETA2 

3 IF ( ABS (THETA3 ) . LT . Cl ) GO TO 5 
I F ( THETA3 . LT . 0 ) GO TO 4 
THETA3 = THETA2-C1 

GO TO 3 

4 THETA3 = THETA2+C1 
GO TO 3 

5 THETA = THETA2 
ENDPRO 

HEAD = THETA *180./C1“ ! 

RUDANG = DELTAC* 180. /Cl 
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FILE: S0M1 



DSL 



A1 



DX 

DY 

X 

Y 



DS 

N 

P 

M 



YERROR 

XERROR 



THETCl 

THETAC 

T 



V*COS ( THETA ) 
V*SIN(THETA) 
INTGRL( 1 . , DX) 
INTGRL ( 0 . , DY) 
TRAJE1-Y 
X-TRAJE2 
XERROR* S IN ( INI) 
DS**2+L**2 
SQRT(N) 

DS/P 
ASIN(M) 
K5*THETC1 
THETAC* 180. /Cl 



PRINT 125 . , RUDANG, HEAD , THECOM, X, TRAJE1 # Y, TRAJE2 , DS 
CONTRL FINTIM=20000 . , DELT=6 . , DELS=20 . 

SAVE ( G1 ) 12 . # HEAD , THECOM 
SAVE ( G2 ) 12 . , TRA JE 1 , Y 
SAVE ( G3 ) 12 . , TRAJE2 , X 
SAVE ( G4 ) 12 . , RUD ANG 

GRAPH (Gl/Gl , DE=TEK618 , P0=0 , .5) TIME ( LE=8 . 0 , SC=2000 . / . 

NI=10 , UN= 1 SEC ' ) , . . . 

HEAD(LE=9,NI=9,L0=10,SC=7. , UN= ' DEGREES ' ) 

THECOM ( LE=9 , NI=9 , P0=8 . , L0=10 , SC=7 . , UN= ' DEGREES ' •) 
GRAPH ( G2/G2 , DE=TEK6 18 , P0=0 , . 5 ) TIME ( LE=8 . , SC=2000 . , . . 

NI=12 , UN= * SEC 1 ) , . . . 

TRAJE1(LE=9,NI=9 / L0=0, SC=6. / UN= f MILES’ ) , . . . 

Y(LE=9 , NI=9 , P0=8 . , L0=0 , SC=6 . ,UN=' MILES' ) 

GRAPH ( G3/G3 , DE=TEK618 , P0=0 , .5) TIME ( LE=8 . , SC=2000 . , . . 

NI=12 , UN= r SEC ' ) 

TRAJE2(LE=9,NI=9,L0=0,SC=8. ,UN=’ MILES 1 ) , . . . 

X( LE=9 , NI=9 , P0=8 . / L0=0 # SC=8. ,UN=’ MILES' ) 

GRAPH (G4/G4 / DE=TEK618 / P0=0 / .5) TIME ( LE=8 . , SC=300 

NI=12 , UN= ' SEC ' ) , . . . 

RUDANG(LE=9,NI=9,LO=-25 / SC=6. , UN= ' DEGREES ' ) 

LABEL (Gl) SHIP HEADING AND DESIRED HEADING VS TIME 

LABEL (Gl) K5= . 75 X=1 IN2=0 . 00349 

LABEL ( G2 ) DESIRED TRAJECTORY AND 

LABEL (G2 ) ACTUAL TRAJECTORY (Y) VS TIME 

LABEL (G3 ) DESIRED TRAJECTORY AND 

LABEL (G3 ) ACTUAL TRAJECTORY (X) VS TIME 

LABEL (G4) RUDDER ANGLE VS TIME 

END 

STOP 
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